home *** CD-ROM | disk | FTP | other *** search
/ Practical Algorithms for Image Analysis / Practical Algorithms for Image Analysis.iso / CH_6.3 / XSGT / KNNS.C < prev    next >
Encoding:
C/C++ Source or Header  |  1999-09-11  |  17.6 KB  |  698 lines
  1. /* 
  2.  * knns.c
  3.  * 
  4.  * Practical Algorithms for Image Analysis
  5.  * 
  6.  * Copyright (c) 1997, 1998, 1999 MLMSoftwareGroup, LLC
  7.  */
  8.  
  9. /*
  10.  * k-N(earest)N(eighbor) S(hells)
  11.  *
  12.  * construct kNN shells on the basis of Voronoi diagram output;
  13.  *
  14.  */
  15. #include <stdio.h>
  16. #include <stdlib.h>
  17. #include "vora.h"
  18.  
  19. #undef    DEBUG
  20. #undef    KNN_DEBUG
  21. #undef MN_DEBUG
  22. #undef    KNN_ARRAY_VER
  23.  
  24. /* globals */
  25. extern int KNN_MAX;
  26. int DISPL_KNNS = 0;             /* when set, display full KNN shells */
  27.  
  28.  
  29. /*
  30.  * for a given Site, evaluate the average, mn, of the nof NN of its nnn NN
  31.  */
  32. float
  33. eval_mn (struct vSite *vS, struct vSite *vsa)
  34. {
  35.   float mn;
  36.   int inn;
  37.  
  38.   mn = (float) 0.0;
  39.   for (inn = 0; inn < vS->nnn; inn++) {
  40. #ifdef MN_DEBUG
  41.     printf ("\t %d (with %d NNs)\n",
  42.             vS->nnsi[inn], (vsa + (vS->nnsi[inn]))->nnn);
  43. #endif
  44.     mn += (float) (vsa + (vS->nnsi[inn]))->nnn;
  45.   }
  46.   mn /= (float) vS->nnn;
  47.  
  48.   return (mn);
  49. }
  50.  
  51.  
  52. /*
  53.  * evaluate, for each site, the average, m, of the nof NN of its NN;
  54.  * ignore sites which have been flagged (eFlag, aoiFlag) out-of-bounds
  55.  *
  56.  * --> Aboav-Weaire law for cellular structures 
  57.  * K.B.Lauritsen, C. Moukarzel & H.J.Herrmann, J. Physique I 3, 1941 (1993)
  58.  */
  59. int
  60. m_of_n (int *ndn, float **m_n, int ns, struct vSite *vsa)
  61. {
  62.   int is, nas;
  63.   int cnn;
  64.   float mn_av;
  65.   struct vSite *vS;
  66.  
  67.   nas = 0;
  68.   for (is = 0; is < ns; is++) {
  69.     vS = (vsa + is);
  70.     if ((vS->eFlag == UnSet) &&
  71.         (vS->aoiFlag == UnSet)) {
  72.       cnn = vS->nnn;
  73.       mn_av = *(m_n[cnn] + ndn[cnn]) = eval_mn (vS, vsa);
  74. #ifdef MN_DEBUG
  75.       printf ("\nM_OF_N: process site: %d\n", vS->sitenbr);
  76.       printf ("\tmn_av: %f, ndn[%d]: %d, *(mn[]+ndn[]): %f\n", mn_av, cnn, ndn[cnn], *(m_n[cnn] + ndn[cnn]));
  77. #endif
  78.       ndn[cnn]++;
  79.       nas++;
  80.     }
  81.   }
  82.   return (nas);
  83. }
  84.  
  85.  
  86. /*
  87.  * evaluate the average area, Ak, of domains in the kNN shells of
  88.  * the current Site
  89.  */
  90. void
  91. eval_aka (struct vSite *vsa, int kNN, double *aka, struct linklist *sha)
  92. {
  93.   int k;
  94.  
  95.   for (k = 0; k < kNN; k++) {   /* loop over k-th NN shells */
  96.     *(aka + k) = xak_kNN_list (vsa, sha + k, k);
  97.   }
  98. }
  99.  
  100.  
  101.  
  102.  
  103. /*
  104.  * evaluate, for the current Site, the product c(0)c(k)/nkNN, as a function
  105.  * of k: c denotes topological charge, n-6, of a given site, nkNN denotes
  106.  * the number of sites in the k-th NN shell; these products represent the
  107.  * correlation function of topological charge as a function of topological
  108.  * distance, k;
  109.  */
  110. int
  111. eval_tctc (struct vSite *vsa, int kNN, double *tctc, struct linklist *sha)
  112. {
  113.   int k;
  114.   double c0, ck_bar;
  115.  
  116.  
  117.   c0 = xc0 (vsa, sha + 0);      /* topol charge of ref Site site */
  118.   *(tctc + 0) = c0 * c0;
  119.   for (k = 1; k < kNN; k++) {   /* loop over k-th NN shells */
  120.     ck_bar = xck_kNN_list (vsa, sha + k, k);
  121.     *(tctc + k) = c0 * ck_bar;
  122.   }
  123.   return ((int) c0);
  124. }
  125.  
  126.  
  127. /*
  128.  * evaluate the fractal measure, c, of a set of kNN shells, associated
  129.  * with a given site, such that N_k = c*k**2, 0<=k<=kNN, N_k representing
  130.  * the nof sites contained within the k-th NN shell, i.e., 
  131.  * N_k = sum from {j=0 to j=k} N_j;
  132.  */
  133. double
  134. eval_frms (int kNN, struct linklist *sha)
  135. {
  136.   int k;
  137.   int nk;
  138.   double frms;
  139.  
  140.   nk = 1;
  141.   for (k = 1; k < kNN; k++) {   /* <= ??? */
  142.     nk += ll_length (sha + k);
  143.   }
  144.   frms = (double) nk / (double) ((kNN - 1) * (kNN - 1));
  145.   return (frms);
  146. }
  147.  
  148.  
  149.  
  150. /*
  151.  * compare list entries 
  152.  *   key: Site index
  153.  */
  154. int
  155. cmpSiteInd (struct nnSite *oldSite, struct nnSite *newSite)
  156. {
  157.   return (SIGN (oldSite->sitenbr - newSite->sitenbr));
  158. }
  159.  
  160. /*
  161.  * init (empty) kNN shell lists
  162.  */
  163. void
  164. init_sk (struct linklist *knnsll)
  165. {
  166.   llzero (knnsll);              /* init empty list */
  167.   llsetsize (sizeof (struct nnSite), knnsll);
  168.  
  169.   llsetmatch (cmpSiteInd, knnsll);
  170. }
  171.  
  172.  
  173.  
  174. /*
  175.  * scan list, comparing index of newSite with that of Sites already 
  176.  * in current list; add Site in order of increasing indices; 
  177.  * employs matching function, here cmpSiteInd();
  178.  * (see also: llsearch() in llist.c )
  179.  */
  180. Boolean
  181. scanSiteInd (struct nnSite *newSite, struct linklist *list)
  182. {
  183.   Boolean MatchStatus = False;
  184.  
  185.   llhead (list);
  186.   do {
  187.     if (((*list->match) (list->clp->item, newSite)) > 0)
  188.       MatchStatus = True;
  189.  
  190.   } while ((MatchStatus != True) && (llnext (list) == True));
  191.  
  192.   return (MatchStatus);
  193. }
  194.  
  195.  
  196. /*
  197.  * insert new site into current list, maintaining the list sorted
  198.  * in order of increasing nnSite indices
  199.  */
  200. void
  201. expand_kNN_list (struct nnSite *nnnS, struct linklist *list)
  202. {
  203.   Boolean MatchStatus;
  204.  
  205.   llhead (list);
  206.   if ((MatchStatus = scanSiteInd (nnnS, list)) == True)
  207.     lladdprev ((char *) nnnS, list);
  208.   else
  209.     lladd ((char *) nnnS, list);
  210. }
  211.  
  212.  
  213. /*
  214.  * remove duplicate entries from kNN list which is assumed to contain
  215.  * kNN site indices in increasing order
  216.  */
  217. void
  218. contract_kNN_list (struct linklist *list)
  219. {
  220.   struct nnSite *nnS, *nnSpp;
  221.  
  222.  
  223.   if (ll_length (list) < 2) {
  224. #ifdef KNN_DEBUG
  225.     printf ("...list contains less than two entries\n");
  226. #endif
  227.   }
  228.   else {
  229.     llhead (list);
  230.     nnS = (struct nnSite *) list->clp->item;
  231.     while (llnext (list) == True) {
  232.       nnSpp = (struct nnSite *) list->clp->item;
  233.       if (nnS->sitenbr == nnSpp->sitenbr) {
  234.         lldelete (list);
  235.         nnS = (struct nnSite *) list->clp->item;
  236.       }
  237.       else {
  238.         nnS = nnSpp;
  239.       }
  240.     }
  241.   }
  242. }
  243.  
  244.  
  245. /*
  246.  * display list
  247.  */
  248. void
  249. show_kNN_list (struct linklist *list, int k)
  250. {
  251.   int is;
  252.  
  253.   printf ("SHOW_KNN_LIST: %2d-NN shell:", k);
  254.   llhead (list);
  255.   if (ll_length (list) != 0) {  /* list not empty */
  256.     is = 0;
  257.     do {
  258.       printf (" %3d ",
  259.               ((struct nnSite *) list->clp->item)->sitenbr);
  260.       is++;
  261.  
  262.       if ((is + 1) % 10 == 0)
  263.         printf ("\n");
  264.     } while (llnext (list) == True);
  265.     printf ("\n");
  266.   }
  267.   else
  268.     printf ("  list empty\n");
  269. }
  270.  
  271.  
  272.  
  273.  
  274. /*
  275.  * extract average area of k-th NN shells from corresponding lists
  276.  */
  277. double
  278. xak_kNN_list (struct vSite *vsa, struct linklist *list, int k)
  279. {
  280.   int is, lk;
  281.   double ak = 0.0;
  282.  
  283.   llhead (list);
  284.   if ((lk = ll_length (list)) != 0) {  /* list not empty */
  285.     do {
  286.       is = ((struct nnSite *) list->clp->item)->sitenbr;
  287.       ak += (double) ((vsa + is)->area);
  288.  
  289.     } while (llnext (list) == True);
  290.     ak /= (double) lk;
  291.   }
  292.   else
  293.     printf ("  list for %d-th NN shell empty\n", k);
  294.  
  295.   return (ak);
  296. }
  297.  
  298.  
  299. /*
  300.  * extract topol charge of ref Site (stored at top of list sha+0)
  301.  */
  302. double
  303. xc0 (struct vSite *vsa, struct linklist *list)
  304. {
  305.   int is, lk;
  306.   double c0 = 100.0;
  307.  
  308.   llhead (list);
  309.   if ((lk = ll_length (list)) != 0) {  /* list not empty */
  310.     is = ((struct nnSite *) list->clp->item)->sitenbr;
  311.     c0 = (double) ((vsa + is)->nnn - 6);
  312.   }
  313.   else
  314.     printf ("  list empty, something wrong\n");
  315.  
  316.   return (c0);
  317. }
  318.  
  319.  
  320. /*
  321.  * extract average topological charge of k-th NN shells 
  322.  * from corresponding lists
  323.  */
  324. double
  325. xck_kNN_list (struct vSite *vsa, struct linklist *list, int k)
  326. {
  327.   int is, lk;
  328.   double tck = 0.0;
  329.  
  330.   llhead (list);
  331.   if ((lk = ll_length (list)) != 0) {  /* list not empty */
  332.     do {
  333.       is = ((struct nnSite *) list->clp->item)->sitenbr;
  334.       tck += (double) ((vsa + is)->nnn - 6);
  335.  
  336.     } while (llnext (list) == True);
  337.     tck /= (double) lk;
  338.   }
  339.   else
  340.     printf ("  list for %d-th NN shell empty\n", k);
  341.  
  342.   return (tck);
  343. }
  344.  
  345.  
  346. /*
  347.  * copy vSite information to root of kNN shell lists
  348.  */
  349. void
  350. init_sk0 (struct vSite *vsa, int sitenbr, struct linklist *list)
  351. {
  352.   struct nnSite nnS, *pnnS = &nnS;
  353.  
  354.   (vsa + sitenbr)->status = InActive;
  355.   pnnS->sitenbr = sitenbr;
  356.  
  357.   expand_kNN_list (pnnS, list);
  358. }
  359.  
  360.  
  361.  
  362. /*
  363.  * set Status of all Sites to Active
  364.  */
  365. void
  366. SetActive (struct vSite *vsa, int ns)
  367. {
  368.   int is;
  369.  
  370.   for (is = 0; is < ns; is++)
  371.     (vsa + is)->status = Active;
  372. }
  373.  
  374.  
  375. /*
  376.  * extract, from Voronoi diagram, kNN shells of a given site, 
  377.  * in order of increasing k: k=0 refers to the Site itself;
  378.  * the NN shells are constructed on the basis of an array of Voronoi sites
  379.  * which gives site coordinates, site index and NN site indices;
  380.  *
  381.  * the desired list of NN shells is implemented as a list of lists, for
  382.  * each site in the set; for each k, there is a list of indices giving the
  383.  * sites of the corresponding shell; for given ROOT site, each site is
  384.  * uniquely assigned to one shell: that is, as k increases, only those
  385.  * site indices area added to the current shell which have not been 
  386.  * assigned to any previous shell; once assigned, a site is given InActive
  387.  * status; the recursion comes to an end when either a maximum k is reached
  388.  * or else, if there are no more Active sites to be added to new shells
  389.  */
  390. struct linklist *
  391. kNNs (struct vSite *vsa, struct linklist *plist, struct linklist *clist, int k)
  392. {
  393.   struct nnSite *pnnS, nnS, *cnnS = &nnS;
  394.   struct vSite *vS;
  395.   int inn;
  396.   Boolean Descend_Status;
  397.  
  398.  
  399. #ifdef DEBUG
  400.   printf ("\nKNNS:\n");
  401.   printf ("...pointers:\n");
  402.   printf ("\tvsa: %p, plist: %p, clist: %p\n", vsa, plist, clist);
  403.   printf ("\textract_kNN: %p\n", kNNs);
  404.   printf ("...new parent list length: %3d\n", ll_length (plist));
  405.   show_kNN_list (plist, k);
  406. #endif
  407.  
  408.   llhead (plist);               /* parent list */
  409.   llhead (clist);               /* child list */
  410.   Descend_Status = False;
  411.   do {                          /* step through parent list */
  412.     pnnS = (struct nnSite *) plist->clp->item;
  413.     vS = vsa + (pnnS->sitenbr);
  414.  
  415.     for (inn = 0; inn < vS->nnn; inn++) {  /* loop over NN ind of vS */
  416.       if ((vsa + (vS->nnsi[inn]))->status == Active) {
  417.  
  418.         cnnS->sitenbr = (vS->nnsi[inn]);
  419.         expand_kNN_list (cnnS, clist);
  420.  
  421.         if (Descend_Status == False)
  422.           Descend_Status = True;
  423.  
  424.         (vsa + (vS->nnsi[inn]))->status = InActive;
  425.       }
  426.     }
  427.   } while (llnext (plist) == True);
  428.  
  429.   if (Descend_Status == True) {
  430.     contract_kNN_list (clist);
  431.     return (clist);
  432.   }
  433.   else {
  434. #ifdef DEBUG
  435.     printf ("\nKNNS: no more Active sites\n");
  436. #endif
  437.     return (clist = (struct linklist *) NULL);
  438.   }
  439. }
  440.  
  441.  
  442.  
  443. /* 
  444.  * code to construct, on the basiis of the Voronoi diagram, the k-NN shells 
  445.  * for each site in a point set 
  446.  *
  447.  * the function construct_kNNs() handles the kNN shells for a single
  448.  * site at a time, and evaluates quantities of interest before going
  449.  * on to the next site; that is, a single array of shells, sk, properly
  450.  * dimensioned to KNN_MAX, is repeatedly used to store the kNN 
  451.  * shells of successively examined sites
  452.  *
  453.  * a note on the data structures:
  454.  *   sk is an array of pointers to struct kNNShell; one such array is
  455.  *   associated with each site, is <= ns, in the set; each entry in sk, 
  456.  *   sk[is], contains the size, kNN, of the actual nof shells built for
  457.  *   a given ``root'' site, and quantities of interest evaluated from that
  458.  *   shell, such as the fractal measure; 
  459.  *   for each root Site, an array, sha, of (maximally) KNN_MAX shells is
  460.  *   constructed (and then recycled for the following root Site): the k-th 
  461.  *   element, sha[k], of that array represents the k-th NN shell of the root 
  462.  *   Site in the form of a list of structs nnSite, containing, among other
  463.  *   entries, the indices (into the Site array vsa) identifying the k-th NN;
  464.  *
  465.  *   if memory consumption still presents a concern, this version could
  466.  *   be modified to dispense with the array, sk, and just store information
  467.  *   for a single Site at a time: in that case, evaluation of quantities 
  468.  *   such as the fractal measure would have to integrated in kNN construction
  469.  *   
  470.  *   a memory-consuming, but conceptually appealing version:
  471.  *   -------------------------------------------------------
  472.  *   for the akNNs version, sk[is] contains, in addition, an array of 
  473.  *   kNN lists; the k-th element, (sk+is)->sha[k], in any such array is 
  474.  *   a list containing the set of indices comprising the k-NN shell of 
  475.  *   the ``root'' or reference site;
  476.  *   sites which have been flagged (eFlag, aoiFlag) out-of-bounds
  477.  *   are not processed: their kNN shells (k>0) contain 0 elements
  478.  */
  479. double
  480. construct_kNNs (int ns, struct vSite *vsa, struct linklist *sha, struct kNNShell *sk)
  481. {
  482. #ifdef KNN_DEBUG
  483.   int j;
  484. #endif
  485.   int is, k;
  486.   int nas;
  487.   double frms;
  488.   struct vSite *vS;
  489.  
  490.   struct linklist *plist, *clist;  /* ptr to parent(k) list */
  491.  
  492.  
  493. /*
  494.  * construct kNN shells, one site at a time
  495.  */
  496.   frms = 0.0;
  497.   nas = 0;
  498.   printf ("\tSite:");
  499.   for (is = 0; is < ns; is++) { /* loop over sites */
  500.     (sk + is)->sitenbr = (vsa + is)->sitenbr;
  501.  
  502.     if (is % 50 == 0)
  503.       printf ("...%3d", is);
  504.     if ((is + 1) == ns)
  505.       printf ("...%3d\n", is);
  506.  
  507. /* set status of all Sites to Active (must be repeated) */
  508.     SetActive (vsa, ns);
  509.  
  510. /* set up empty list structs */
  511.     for (k = 0; k < KNN_MAX; k++)
  512.       init_sk (sha + k);
  513.  
  514. /* loop over successive k-NN shells of current ``root'' site */
  515.     init_sk0 (vsa, (vsa + is)->sitenbr, sha + 0);
  516.     if (DISPL_KNNS == 1)
  517.       show_kNN_list (sha + 0, 0);
  518.  
  519.     if (((vS = (vsa + is))->eFlag == UnSet) &&
  520.         (vS->aoiFlag == UnSet)) {  /* Site out-of-bounds ? */
  521.       k = 0;
  522.       plist = (struct linklist *) (sha + 0);
  523.       do {
  524.         if (++k >= KNN_MAX) {
  525. #ifdef KNN_DEBUG
  526.           printf ("...reached max NN shell: ");
  527.           printf (" k: %d\n", k);
  528. #endif
  529.           plist = (struct linklist *) NULL;
  530.         }
  531.         else {
  532.           clist = (struct linklist *) (&sha[k]);
  533.           plist = kNNs (vsa, plist, clist, k);
  534.  
  535. #ifdef KNN_DEBUG
  536.           printf ("...%d-th shell completed\n", k);
  537.           for (j = 0; j <= k; j++)
  538.             show_kNN_list (sha + j, j);
  539. #endif
  540.         }
  541.       } while (plist != (struct linklist *) NULL);
  542.       (sk + is)->kNN = k;
  543.  
  544.       if (DISPL_KNNS == 1) {
  545.         printf ("...%d-NN shells for Site %d\n",
  546.                 k, vS->sitenbr);
  547.         for (k = 0; k < (sk + is)->kNN; k++)
  548.           show_kNN_list (sha + k, k);
  549.       }
  550.  
  551. /*
  552.  * employ kNN shells of current root Site to evaluate further
  553.  * quantities of interest; soter these in sk[is]
  554.  *
  555.  * fractal measure:
  556.  *   c, in Nk = c*k**2, where Nk denotes the number 
  557.  *   of domains in the k-th NN shell of a given site
  558.  */
  559.       (sk + is)->frms = eval_frms ((sk + is)->kNN, sha);
  560. #ifdef KNN_DEBUG
  561.       printf ("\n...fractal measure (from %d shells): %lf\n",
  562.               (sk + is)->kNN, (sk + is)->frms);
  563. #endif
  564.       frms += (sk + is)->frms;
  565.       nas++;
  566.  
  567.  
  568. /*
  569.  * average area, <Ak>, of domains in k-th NN shell, as function of k
  570.  */
  571.       eval_aka (vsa, (sk + is)->kNN, (sk + is)->aka, sha);
  572.  
  573. /*
  574.  * correlation function, tccf(k) = c(0)c(k), of topological charge,
  575.  * c=n-6, as function of topological distance, k, for given site;
  576.  */
  577.       (sk + is)->tc = eval_tctc (vsa, (sk + is)->kNN, (sk + is)->tctc, sha);
  578.     }
  579.     else {                      /* Site skipped */
  580.       (sk + is)->kNN = 0;
  581.     }
  582.   }
  583.   frms /= (double) nas;
  584.   return (frms);
  585. }
  586.  
  587.  
  588.  
  589. /* 
  590.  * code to construct k-NN shells for each site in a point set 
  591.  * on the basis of Voronoi diagram:
  592.  * this version, construct_akNNs() handles an entire array of shells,
  593.  * one for each site in the set (see loop over sites); this has the 
  594.  * disadvantage of requiring large amounts of memory, but has the
  595.  * advantage of making available all requisite data for subsequent
  596.  * computation of topological quantities
  597.  *
  598.  * 3/94: not completely debugged, not maintained
  599.  *
  600.  * a note on the data structures:
  601.  *    sk is an array of pointers to struct kNNShell; one such array is
  602.  *    associated with each site, is <= ns, in the set; each entry in sk, 
  603.  *    sk[is], contains the size, kNN, of the actual nof shells built for
  604.  *    a given ``root'' site, and an array of kNN lists; the k-th element, 
  605.  *    (sk+is)->sha[k], in any such array is a list containing the set 
  606.  *    of indices comprising the k-NN shell of the ``root'' site
  607.  */
  608. #ifdef KNN_ARRAY_VER
  609. struct kNNShell *
  610. construct_akNNs (int ns, struct vSite *vsa)
  611. {
  612.   int is, j, k, kNN;
  613.  
  614.   struct kNNShell *sk;          /* array of struct kNNShell */
  615.   struct linklist *plist, *clist;  /* ptr to parent(k) list */
  616.   struct linklist *sksi;
  617.  
  618. /*
  619.  * mem alloc
  620.  */
  621.   if ((sk = (struct kNNShell *) calloc (ns, sizeof (struct kNNShell))) == NULL)
  622.       fail_alloc ("sk", 1);
  623.  
  624.   for (is = 0; is < ns; is++) { /* loop over sites */
  625.  
  626. /*
  627.  * handle shell structure mem alloc
  628.  */
  629.     if (((sk + is)->sha = (struct linklist *) calloc (KNN_MAX, sizeof (struct linklist))) == NULL)
  630.         fail_alloc ("(sk+is)->sha", 1);
  631.  
  632.     sksi = (sk + is)->sha;      /* array of kNN lists */
  633.  
  634. /*
  635.  * printf("\n...size of struct linklist: %d\n", 
  636.  * sizeof(struct linklist));
  637.  * printf("...size of sksi: %d\n", sizeof(sksi));
  638.  */
  639.  
  640. /*
  641.  * set status of all Sites to Active (must be repeated)
  642.  */
  643.     SetActive (vsa, ns);
  644.  
  645. /* set up empty list structs */
  646.     printf ("\n...initialize empty k-NN shell lists\n");
  647.     for (k = 0; k < KNN_MAX; k++)
  648.       init_sk (&((sk + is)->sha[k]));
  649.  
  650. /* loop over successive k-NN shells of current ``root'' site */
  651.     init_sk0 (vsa, (vsa + is)->sitenbr, &((sk + is)->sha[0]));
  652.     show_kNN_list (&((sk + is)->sha[0]), 0);
  653.  
  654.     k = 0;
  655.     plist = &((sk + is)->sha[0]);
  656.     do {
  657. /*                      clist = (struct linklist *)(sksi+k); */
  658. /*                      
  659.  * printf("\n...pointers:\n");
  660.  * printf("\tsksi: %p, sksi+k: %p, ++sksi: %p\n",
  661.  * sksi, sksi+k, ++sksi);
  662.  * printf("\tplist: %p, clist: %p\n", plist, clist);
  663.  * --sksi;
  664.  */
  665.  
  666.       if (++k >= KNN_MAX) {
  667.         printf ("...reached max NN shell, k: %d\n", k);
  668.         plist = (struct linklist *) NULL;
  669.       }
  670.       else {
  671.         printf ("...new k: %d\n", k);
  672.         clist = &((sk + is)->sha[k]);
  673.  
  674.  
  675.         plist = kNNs (vsa, plist, clist, k);
  676.  
  677. #ifdef KNN_DEBUG
  678.         printf ("...%d-th shell completed\n", k);
  679.         for (j = 0; j <= k; j++)
  680.           show_kNN_list (&((sk + is)->sha[j]), j);
  681. #endif
  682.       }
  683.     } while (plist != (struct linklist *) NULL);
  684.     (sk + is)->kNN = kNN = k;
  685.  
  686.     printf ("...%d-NN shells for Site %d\n", k, (vsa + is)->sitenbr);
  687.     for (k = 0; k < (sk + is)->kNN; k++)
  688.       show_kNN_list (&((sk + is)->sha[k]), k);
  689.  
  690.     if (((sk + is)->sha = (struct linklist *) realloc ((sk + is)->sha, kNN * sizeof (struct linklist))) == NULL)
  691.         fail_alloc ("realloc (sk+is)->sha", 1);
  692.   }
  693.   for (is = 0; is < ns; is++)
  694.     free ((sk + is));
  695.   return (sk);
  696. }
  697. #endif
  698.